Non-linear screening of external charge by doped graphene 
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We solve a nonlinear integral equation for the electrostatic potential in doped graphene due to an external 
charge, arising from a Thomas-Fermi (TF) model for screening by graphene's n electron bands. In particular, 
we study the effects of a finite equilibrium charge carrier density in graphene, non-zero temperature, non-zero 
gap between graphene and a dielectric substrate, as well as the nonlinearity in the band density of states. Effects 
of the exchange and correlation interactions are also briefly discussed for undoped graphene at zero temperature. 
Nonlinear results are compared with both the linearized TF model and the dielectric screening model within 
random phase approximation (RPA). In addition, image potential of the external charge is evaluated from the 
solution of the nonlinear integral equation and compared to the results of linear models. We have found generally 
good agreement between the results of the nonlinear TF model and the RPA model in doped graphene, apart 
from Friedel oscillations in the latter model. However, relatively strong nonlinear effects are found in the TF 
model to persist even at high doping densities and large distances of the external charge. 

PACS numbers: 73.22.Pr 

Keywords: graphene, Thomas-Fermi model, screening, charged impurity, Friedel oscillations 



I. INTRODUCTION 

Over the period of just five years since its first incep- 
tion in the laboratory, 1 graphene has developed into one of 
the currently most active research areas in the nano-scale 
physics. 2 One of the most important, and certainly most elu- 
sive, problems in graphene research is concerned with its elec- 
trical conductivity, especially in the regime close to zero dop- 
ing of graphene, where the conductivity exhibits a peculiar 
minimum. 3 ' 4 - 5,6,7 Besides several other scattering mechanisms 
for charge carriers in graphene, it is believed that a special 
role in graphene's conductivity is played by the carrier scatter- 
ing on charged impurities, which are ubiquitous in graphene's 
surroundings. In that context, significant progress has been 
achieved in understanding the conductivity of graphene by us- 
ing the Boltzmann transport theory for charge carrier scatter- 
ing on linearly screened charged impurities within the random 
phase approximation (RPA) i 8 ' 9 ' 10 However, because of the re- 
duced dimensionality, and especially because of the semi- 
metallic nature of graphene's ir electron bands, the problem 
of screening of charged impurities remains open. In that con- 
text, other approaches have also been undertaken, including 
a full scattering theoretical treatment of Coulomb impurities 
embedded within the graphene plane ^ 11 ' 12 : 13 ' 14 as well as non- 
linear screening of external charges studied by means of the 
Thomas-Fermi (TF) ; 15 ' 16,17 ' 18 Thomas-Fermi-Dirac (TFD), 19 
and Density Functional Theoretical (DFT) schemes. 20 

While graphene's applications in nanoelectronics are pri- 
marily concerned with charged impurities trapped in an in- 
sulating substrate^LS 2 . screening of external charges is also 
of interest for sensor applications of graphene in detecting 
atoms or molecules^ which may be either adsorbed on the 
upper surface of graphene^i^ 5 . or intercalated in the gap be- 
tween the graphene and the substrate. 26 Further applications 
include image-potential states of electrons near graphene* 2 - 7 ^ 
as well as the image and friction forces on slowly moving 
ions that may affect the kinetics of chemical reactions tak- 
ing place in the vicinity of grapheneji^ 2 ^ All these aspects 



of screening of external charges by graphene are expected to 
be strongly influenced by the presence of nearby dielectric 

materials, 30 . 31 . 32 ^ 3 . 34 

One of the most important issues in theoretical studies of 
screening of external charges is concerned with applicabil- 
ity of the linear response theory for intrinsic, or undoped 
graphene. Namely, with its valence and conducting tt elec- 
tron bands touching each other only at the K and K' points 
of the Brillouin zone, 2 graphene behaves as a zero-gap semi- 
conductor, so that its polarizability is greatly reduced when its 
Fermi level lies close to the neutrality (or Dirac) point char- 
acterizing the regime of zero doping. In that context, it was 
shown within the RPA approach that screening of external 
charges by intrinsic graphene at zero temperature is charac- 
terized merely by a renormalization of graphene's background 
dielectric constant due to inter-band electron transitions i 8,35,36 
However, when graphene is doped up to a certain number den- 
sity n (per unit area) of charge carriers, e.g., by applying an 
external gate potential, then its Fermi level shifts away from 
the neutrality point and the linear screening is expected to be- 
come appropriate, even at zero temperature. It is therefore 
desirable to determine the parameter range where nonlinear 
effects in screening of an external charge set in by contrast- 
ing the results from linear screening models with those from 
available nonlinear models, such as TF and DFT. 

In that context, Katsnelson 16 and Fogler et alJl have solved 
the nonlinear TF model, first proposed by DiVicenzo and 
Melsi£ for intrinsic graphene (i.e., n = 0) in the presence 
of an external point charge. These authors found unusually 
long ranged induced density of charge carriers in the plane 
of graphene, 17 and showed that the linear approximation to 
the TF model for the induced potential is likely to overesti- 
mate the contribution of scattering on charged impurities to 
the resistivity of graphene. 16 However, performance of the TF 
model has been recently criticized for intrinsic graphene in 
the presence of sufficiently weak periodic perturbations vali- 
dating linear screening within the RPA. 37 On the other hand, 
the above nonlinear TF model, augmented by the exchange (or 
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Dirac) interaction in the local density approximation (LDA), 
proved to be valuable in estimating the effective potential fluc- 
tuations in doped graphene due to randomly distributed mul- 
tiple charged impurities^ A similar problem in the presence 
of multiple charged impurities was also tackled by a more ad- 
vanced DFT approach including both the exchange and cor- 
relation (XC) interactions in LDA. 20 All the above models 
were formulated assuming a zero temperature, linear density 
of states (DOS) of the ir electron bands, and a zero gap be- 
tween graphene and substrate. 



In this paper, we take up the simple TF model for a single 
point charge Ze, a distance zq away from graphene , 15 ' 16 ' 17 ' 37 
and generalize it to include the effects of a non-zero ground- 
state charge carrier density n, a non-zero temperature T, 
and the presence of a substrate at a non-zero distance h 
from graphene.— We assume that the external charge is 
weak/distant enough to have negligible effects on the struc- 
tures of graphene's DOS, apart from its shift due to local 
charging of graphene, but we allow for large displacements of 
the Fermi level away from the neutrality point by including the 
nonlinear corrections to the DOS in our model 2 . By varying 
the magnitude \n\, we are able to examine the effects of dop- 
ing, whereas any dependence on the sign of n will be a signa- 
ture of nonlinear effects in screening by graphene. (Note that 
changing the sign of n with fixed sign of the external charge Z 
in the TF model is equivalent to changing the sign of Z with 
fixed sign of n.) 



We perform a series of numerical solutions of the nonlinear 
integral equation resulting from the TF model for the in-plane 
value of total electrostatic potential for a range of values of 
n and zq, for both zero and room temperature, in the cases 
of both free graphene and an Si02 substrate with the gaps 
h = and 1 A. In a special case of free, intrinsic graphene at 
zero temperature, we also solve the nonlinear TF model aug- 
mented by the XC energy terms of Polini et al^ in order to 
estimate the importance of the exchange and correlation in- 
teractions within the TF approach to screening of an external 
charge. While the results obtained for the radial dependence 
of the in-plane potential could be directly used to discuss non- 
linear effects in graphene's conductivity within the Boltzmann 
transport theory, we turn our attention in the present work to 
using our numerical solutions of the TF model to evaluate the 
nonlinear image potential of an external charge, which pro- 
vides an integrated measure of graphene's screening ability 
and is also of interest in recent studies of the electron image 
states. 27,28 Finally, we compare our nonlinear results for both 
the in-plane potential and the image potential with those from 
the linearized TF (LTF) model and the temperature dependent 
RPA dielectric-function approach. 



After outlining the basic theory in section 2, we discuss 
our results in section 3, and present our concluding remarks 
in section 4. Note that gaussian electrostatic units are used 
throughout unless otherwise explicitly indicted. 



H. THEORY 

We use a cartesian coordinate system with coordinates 
{r, z}, where r = {x, y}, and assume that graphene is placed 
in the z — plane. A semi-infinite substrate with dielectric 
constant e s is assumed to occupy the region z < —h under- 
neath the graphene, whereas the region z > —his assumed to 
be vacuum or air. 18 We assume that the ground state of such a 
system, under the gating conditions at temperature T, is char- 
acterized by a uniform density per unit area of charge carriers 
in the graphene, given by 



n(p) = / de p(e) 
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where p(e) is the DOS in graphene's ir electron bands given 
in Ref.— , = (/c^T) - , and p is the chemical potential. Note 
that for electron (hole) doping, one has n > (n < 0) and 
consequently p > (p < 0), whereas intrinsic graphene is 
characterized by n = and p = 0. For sufficiently low 
doping levels, such that, e.g., \p\ < 1 eV, one may use the 
linearized band DOS 2 ^ ~ — ^ 
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where gd = 4 is the spin and valley degeneracy factor, v f the 
Fermi speed of graphene which we take to be w c/300, where 
c is the speed of light in vacuum, and dilog is the dilogarithm 
function.— 

We wish to evaluate the total electrostatic potential in the 
system, $(r, z), due to an external point charge Ze placed at 
a fixed position {0, zo}, where e is the charge of a proton. 
This perturbation will induce surface charges on the surface 
of the substrate and on the graphene with the densities per 
unit area cr S ub(r) and cr gr (r), respectively. 18 Using the tilde to 
denote the Fourier transform with respect to coordinates in the 
graphene plane, r — > k, we can write the total potential as the 
sum $ = <l> GX t + ^ind, where 
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is the potential of the external charge screened by the di- 
electric constant, e^, of the "host" environment in which that 
charge resides (eh = lforzo > —h and eh = e s forzo < — h), 
and 
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is the total induced potential in the system. Next, one can 
eliminate the Fourier transform of the charge density on the 
substrate, <7 su b(k), by using the boundary condition^ 2 - 
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dz 
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and write for the total induced potential 



$ind(k, z) = — cr gr (k) I e 
2n Ze e s — 1 
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where sgn is the signum function. 



sgn(z + /iX6) 



A. Nonlinear TF model 

In the spirit of a temperature-dependent TF model, we ex- 
press the induced charge density in graphene a a 18 ' 43 ' 44 



and hence — 1. Note that the integral equation, Eq. (0 
with Eq. ( [Tol l, implies an asymmetry with respect to the sign 
of Zo when h > 0, which is lost in the zero gap case. 

In order to discuss the effects of XC interactions within 
the TF model, we note that density dependent expressions 
for both the exchange and correlation energy per electron in 
graphene, are available in the LDA only for density variations 
with respect to the equilibrium case of intrinsic, or undoped 
graphene having /x = 0, in the limits of zero temperature, 
zero gap, and linearized band DOS. 19,20 Therefore, we spe- 
cialize Eq. (0 to those parameters and convert it to an integral 
equation for the potential energy U (r) = e<j>(r), 
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where n (fi) is given by Eq. (Q]), and where we have denoted 
the total electrostatic potential in the graphene plane by 



(r) = $(t,z) 



(8) 



By using the inverse Fourier transform of Eq. © in the ex- 
pression $ = $ cxt + <E>ind in which we set z = 0, we obtain 
the following nonlinear integral equation for </>(r), 18 
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where 
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is the value of the potential due to the external charge in the 
presence of substrate alone, evaluated at z — 0. 

We further convert Eq. (0 with Eq. ( [Tol l into an integral 
equation for the potential energy, defined by U(r) — e</>(r), 
and solve it numerically for a range of the model parame- 
ters, as discussed in the following section. Owing to the axial 
symmetry of the problem, an angular integral can be read- 
ily completed in Eq. (0, giving a one-dimensional integral 
equation for U(r), which is particularly difficult to solve be- 
cause of the singular nature of its integrand, especially for 
intrinsic graphene at zero temperature.^*^ We have mapped 
the interval r e [0, oo) onto a finite interval, partitioned the 
function U at up to 2400 (typically 800) points, and used the 
f solve routine while regularizing the integrand. As a check 
of our method for free graphene, we substituted the solution 
U(r) = e<p(r) into Eq. (0 and verified that its spatial integral 
yields — Ze. 

We note that a more compact form of the integral equation, 
Eq. © with Eq. ( [Tol l, may be obtained for a zero gap (h = 0) 
between graphene and the substrate, giving rise to an overall 
effective background dielectric constant e° g = (e s + 1) /2, as 
is usually done in the literature on graphene i 16 ' 19 i 20 i 35 ' 36 In that 
case, the free graphene limit is recovered by setting e s = 1 
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where K is the complete elliptic integral of the first kind, 40 
and the function n(U) is obtained from Eq. (ffj) in the limit of 
zero temperature as follows 
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For the XC potential energy, V xc (n), we use the expressions 
provided by Polini et al. 20 in their subsections II. A and B. As 
a reference, we quote here the dominant contribution to this 
energy at low densities, which is on the order of 



Kc(n) ^ — - 



\n\ ln( r-^r 

\n\ 
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10) for |n| <^ no = rj x 0.7635 A~ 2 , where 77 is a cut-off pa- 
rameter of the theory taking a value from the interval (0, l]. 20 
Finally, we note that Eq. (fTTT i with V xc (n) set to zero appeared 
in previous studies using the TF model . 15 i 16 ' 17 ' 37 



B. Linear models 

Going back to the TF integral equation, Eq. (0, if the to- 
tal potential </>(r) in the plane of graphene may be treated as 
a weak perturbation of the equilibrium carrier charge density, 
then one can linearize Eq. ©, <r gr (r) rj — e 2 0(r)<9n(/x)/9/i, 
with n(fi) defined in Eq. ([TJ. This facilitates the use of the 
Fourier transform in solving Eq. (0 thereby giving an approx- 
imate expression for the total potential 
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is the background dielectric constant due to substrate, 
v c(k) — 27re 2 /fc, and LT(fc) is the polarization function of 
free graphene, which is constant in the LTF model, given by 
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IItf = dn(n)/dfj,. In Eq. (fl4b . one needs to use the Fourier 
transform of the potential in Eq. ( fTOb , which is given by 



0o (k) 
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-h< z < 
zq < — /;, 



where e° g = e bg (0) = (e s + l)/2. 

In the zero gap limit, one obtains from Eq. ( TBI a more 
compact expression for the total potential in the LTF model, 
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It is clear then that, at zero temperature, intrinsic graphene 
cannot screen external charges in the LTF model because 
q s — > 0. 8 On the other hand, when either n ^ or T > 0, 
the inverse Fourier transform of Eq. ( fTTb gives a total poten- 

(16jjial with the asymptotic formiS 4>(r) ~ (^Zee^^j / (q 2 r 3 ) 

for r ^> qj 1 ^ |j»o|, an d with the limiting value at the 

origin 0(0) = [Ze/ (e£ g *o)] [l - Ce c E x (C)], where C = 

9s z o/ e b<T ar, d Ei is the exponential integral function^ 



(17) 



where the inverse screening length of free graphene, q s — 
27re 2 IlTF, is obtained from Eq. (O within the linearized DOS 

as£ 



In [2 cosh(/3/x/2)] 



The expression (TPfl) can also be used to obtain the total po- 
tential based on the RPA model if one employs the polariza- 
tion function, which is obtained for the linearized DOS at non- 

„8.38.39 



(18) zero temperature as 
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where we have defined a thermal inverse screening length by 
q t = 2/(PHvp). Note that /i, which is to be used in Eq. ( fT9l ), 
can be obtained from Eq. for any given temperature and 
equilibrium charge carrier density n. In the zero temperature 
limit, [i — > £p where ep = hvpkpsgn(n) is the Fermi energy 
with kp = y/Tr\n\ being the Fermi momentum in graphene 
with the equilibrium charge carrier density n, so that one ob- 
tains from Eq. (Jjpa^ 
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where H is the Heaviside unit step function. Unlike the LTF 
case, we see that IIrpa^) = kj (AHvf) in intrinsic graphene 
at zero temperature. Since this is also the short wavelength 
limit of IIrpa(A:) when n ^ 0, one may assert that the 
RPA result will yield a value for the total potential that is 

reduced by an approximate factor of 1 - 



2f° 



where r s = e / (hvp) s» 2.2, when compared to the corre- 
sponding value from the LTF approach for kp^/r 2 + Zq <C 1 
at zero temperature and zero gap. On the other hand, one can 
expect that the total potential will exhibit Friedel oscillations 
for kpT 3> 1 due to non-analyticity of the RPA polarization 
function ( f20T > at k = 2kp, which will be gradually dampened 
as the temperature increases^ 



C. Image interaction 

Once the integral equation, Eq. (0, is solved for the total 
potential in the plane of graphene, one can use Eq. (0 to eval- 
uate the induced charge density in graphene, whose Fourier 
transform may be used in Eq. (O to yield the total induced 
potential for any value of z. This may be then used to calcu- 
late the nonlinear image force on the external charge from the 
definition 
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Once the zq dependence of the image force is determined, the 
corresponding image potential may be obtained from the defi- 
nition Vi m (zo) = J z °° dz' Fi m {z' ). While in the nonlinear TF 
case this integration has to be performed numerically, in a lin- 
ear theory one may use instead the usual definition of image 
potential as a classical self-energy^ VS m (z ) = |Ze$ ind (r = 
0, z = Zq), which gives for z$ > 0^ 

\2 



V im (z ) = 
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By using the LTF model, where vcik)!^^ 
obtains in the zero gap case 

2 [l- £ « g -2Ce 2 «E 1 (2C)], 
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is worthwhile men- 
gives asymptotically 
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V im ~ - (Ze) 2 [1/ (4z ) - 1/ (8q s z 2 )] for a heav- 
ily doped graphene and/or sufficiently large distance, 
such that q s ZQ ^> 1, On the other hand, in the op- 
posite limit, q s ZQ -C 1, one finds to the leading order 

V im ~ {Ze) 2 fl/egg - l) / (4s ), as if graphene were 
totally absent. When the RPA polarization function at zero 
temperature, Eq. d20l i, is used in (f22b in the zero gap case, one 
can show that similar limiting forms of the image potential 
exist, except that the effective background dielectric constant, 
e£ g , is to be replaced by e£ g + irr s /2 w e" g + 3.44 when 

kpz < 1^ 

III. RESULTS 

We first consider the case of a positive charge with Z = 1 
a distance 2 A above graphene lying on an SiC>2 substrate 



(e B = 3.9) with several gap heights h, several equilibrium 
charge carrier densities n, and zero temperature. This situa- 
tion may be representative of a Li atom adsorbed atop sup- 
ported graphene, where the effective charge transfer is found 
to be around Z — 0.9, whereas the local DOS exhibits a res- 
onant feature at about 0.9 eV above the neutrality point of 
graphene's it electron band due to hybridization with lithium's 
2s orbital. 24 Besides undoped graphene with n = 0, which 
was studied previously ; 15 ' 16 ! 17 we also analyze the cases of 
both electron (n > 0) and hole (n < 0) doping of graphene 
by a gate potential, making sure that the Fermi level stays 
well below any chemisorption resonances in graphene's DOS 
(n < 10 13 cm -2 for Li atom^) 



In Fig. 1 we show in the left column (1) the results for the 
potential energy U(r) — ecf>(r), with <j>(r) obtained from the 
nonlinear TF equation Eq. (O at zero temperature for n = 
(the upper thick solid line), ±10 12 (thin solid and dashed 
lines, respectively), and ±10 13 cm~ 2 (the lower thick solid 
and thick dashed lines, respectively), and with h = (pan- 
els a), 1 A (panels b), and oo (i.e., free graphene, panels c). 
For the purpose of comparison, we also show in the right col- 
umn (2) of Fig. 1 the corresponding results obtained from 
both the LTF (dash-dotted lines) and the RPA (dotted lines) 
models (with the line thicknesses matching those in the left 
column), with (f>(r) calculated from Eq. ( TT4l using the ap- 
propriate polarization functions at zero temperature. [As a 
reference, note that, for free graphene, the LTF result with 
n = actually shows the value of the unscreened potential 
in the plane of graphene, Uo(r) = e</>o(r) with 4>q(t) given 
in Eq. | [T0t . whereas the corresponding RPA result shows that 
same potential reduced by the dielectric constant of intrinsic 
graphene, 1 + 7rr s /2 w 4.44.] One can see in Fig. 1 that the 
main effects on the potential come from increasing the doping 
density \n\. While all models exhibit strong variation with n at 
large distances r, one notices that both the nonlinear TF and 
the RPA results are surprisingly concentrated in a relatively 
narrow range of values for the potential at short distances for 
all densities n. This seems to corroborate conclusions from 
a DFT study that the induced density variations in graphene 
seem to saturate with increasing level of doping. 20 

While the LTF model appears to be a rather poor approx- 
imation to the nonlinear TF results at short distances r, their 
agreement improves at large distances with increasing density 
|n|, as expected. Most strikingly, the RPA model gives a sur- 



prisingly good approximation to the nonlinear TF results at 
short distances for all densities n, while exhibiting Friedel os- 
cillations around the LTF results at large distances for n ^ 0, 
with wavelengths that obviously scale with However, 
for n = 0, one sees an increasing disagreement between the 
nonlinear TF and the RPA models with increasing distance, 
which may be attributed to a poor performance of the TF 
model in intrinsic graphene for induced charge carrier den- 
sities below 10 11 cm~ 2 , as suggested recently by Brey and 
Fertig. 37 On the other hand, the TF model presumably gives a 
correct order of magnitude for nonlinear effects, if any, when 
the doping density increases, which are best seen by ana- 
lyzing the effect of changing the sign of n (equivalently, the 
sign of Z), because linear models are insensitive to this sign. 
In that respect, one can clearly notice in the left column of 
Fig. 1 differences between the potentials U+(r) for n > and 
U- (r) for n < in the nonlinear TF model, which will be 
further discussed in Fig. 3 below. 

Finally, one notices in Fig. 1 that, while the presence of 
a finite gap between graphene and substrate does not affect 
qualitative behavior of the results, its quantitative effects may 
not be neglected in the values of the potential for all densities 
shown. While this is particularly obvious at short distances 
for the nonlinear TF results, it is also interesting to see how 
Friedel oscillations in the RPA model increase in amplitude 
with increasing gap. In fact, we have found that the RPA po- 
tential can even change its sign at large distances r for free 
graphene with large enough \n\. Given that the size of gap 
is poorly defined parameter, with a plausible value of around 
h = 1 one should be aware of its role in the total po- 

tential in graphene due to external charges. 



We next consider in Fig. 2 graphene on an Si02 substrate (T = 300 K, panels c and d) temperatures, with a charge 
with the gap h = 1 A, both at zero (panels a and b) and room 




FIG. 1 : The potential energy, U (r) = e<f>(r) (in eV), due to an external proton at distance zo = 2 A above graphene at zero temperature, as a 
function of the radial distance r (in A) in the plane of graphene lying on an Si02 substrate with the gap heights h = (panels a), 1 A (panels 
b), and oo (free graphene, panels c). Results from the nonlinear TF model are shown in the column 1 for equilibrium densities n — (upper 
thick [black] solid line), ±10 12 (thin [red] solid and dashed lines, respectively), and ±10 13 cm" 2 (lower thick [blue] solid and dashed lines, 
respectively). Results from the linearized TF model and the RPA model are shown, respectively, by dash-dotted and dotted lines in the column 
2 for densities \n\ = (upper thick [black] lines), 10 12 (thin [red] lines), and 10 13 cm -2 (lower thick [blue] lines). 



Z = 1 placed at larger distances of z$ = ±10 A away from 
graphene. With zq = 10 A (panels a and c) we can represent a 
distant charge above graphene, such as a slowly moving ion, 29 
or an electron in an image-potential state, 28 whereas the case 
zo = —10 A (panels b and d) represents a technologically 
relevant case of a charged impurity trapped deep in the SiC>2 
substratei^ 2 . We compare the nonlinear TF results with those 



from the RPA model for |n| = 0, 10 12 , and 10 13 cm" 2 , shown 
with the same line styles and thicknesses as in Fig. 1. While 
the RPA results seem to be quite close, apart from the Friedel 
oscillations, to those of the nonlinear TF model for n > 0, 
the agreement between those two models seems to have wors- 
ened at short distances for n = when compared to Fig. 1, 
which may have to do with the problematic performance of 
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' lower: 



r[A] 



FIG. 2: The potential energy, U (r) = e(f>(r) (in eV), due to an external proton at distances zo = ±10 A (left and right columns, respectively) 
from graphene at T — (top row) and T = 300 K (bottom row), as a function of the radial distance r (in A) in the plane of graphene lying 
on an Si02 substrate with the gap height h = 1 A. Results from the nonlinear TF model are shown for equilibrium densities n — (upper 
thick [black] solid line), ±10 12 (thin [red] solid and dashed lines, respectively), and ±10 13 cm -2 (lower thick [blue] solid and dashed lines, 



respectively). Results from the RPA model are shown by dotted lines for densities 
and 10 13 cm -2 (lower thick [blue] line). 



(upper thick [black] line), 10 (thin [red] line), 



the nonlinear TF model in intrinsic graphene exposed to weak 
perturbations, as mentioned previously. 37 

On the other hand, one notices in Fig. 2 a much greater 
spread in the relative magnitudes of the potential at short dis- 
tances than in Fig. 1. This is partly due to the effect of doping 
in the presence of a much weaker external perturbation in Fig. 
2 than in Fig. 1, so that the induced density variations involved 
in the results in Fig. 2 have not reached the effect of satura- 
tion mentioned earlier. 20 An another cause for a larger spread 
of the potential at short distances in Fig. 2 comes from the 
nonlinear effects, which will be further discussed in Fig. 3. 



As regards the effect of finite temperature, one notices that 
its main role is to dampen the potential in intrinsic graphene 
at distances r > 10 A, both in the nonlinear TF and the RPA 
cases. This can be explained by assessing the TF inverse 
screening length in Eq. ([T8l > in the zero density and the zero 
temperature limits, giving q s — > Ar s q t In 2 and q s — > 4r s kp, 
respectively. Therefore, one may conclude that screening of 
the potential at large distances due to a non-zero temperature 
will prevail only for low enough charge-carrier densities, such 
that |n| < [21n2fc B T/ {hv F )f /it w 10 11 cm" 2 atroom tem- 
peratures. [In fact, we have checked that nonlinear TF results 
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for |n| = 10 11 cm~ 2 at zero temperature (not shown) are quite 
close to those shown in Fig. 2 for intrinsic graphene at room 
temperature.] The effects of temperature on the nonlinearity 
of the potential are further discussed in Fig. 3. On the other 
hand, while the Friedel oscillations are still visible in Fig. 2 in 
the RPA results for zero temperature at large distances r for 
n 7^ 0, although they seem to be reduced in relative amplitude 
by the increased distance |zo| when compared to the oscilla- 
tions seen in Fig. 1, one notices that the increased temperature 
dampens the Friedel oscillations in Fig. 2, as expected. 

We note finally that, by analyzing the asymmetry in the re- 
sults with respect to the change in sign of zq in Fig. 2, we yet 
again emphasize the role of a finite gap, because all results 
would be independent of that sign in the zero gap case. It is 
remarkable that the gap of only h = 1 A affects not only the 
values of the potential ar short distances, but also the magni- 
tudes of the asymmetry in the nonlinear TF results with re- 
spect to the sign of n ^ at short distances. 

Nonlinear effects in screening of an external charge by 
doped graphene, seen in Figs. l(bl) and 2(a), are summarized 
in Fig. 3, with the inclusion of the results for doping density 
of \n\ = 10 11 cm~ 2 . We show the ratio U-{r)/U + {r) of the 
potential energies [/_ (r) and U+(r), which are obtained from 
Eq. (0 with, respectively, negative (hole doping) and positive 
(electron doping) signs of densities |n| = 10 11 (solid lines), 
10 12 (dashed lines), and 10 13 cm -2 (dash-dotted lines), for a 
charge Z = 1 at two distances with two temperatures: zo = 2 
A and T = (panel a), z = 10 A and T = (panel b), 
and zo = 10 A and T = 300 K (panel c), for graphene lying 
on an SiC>2 substrate with the gap h = 1 A. One notices in 
Fig. 3 that the ratio U- (r)/U+(r) may reach quite large val- 
ues (up to a factor of two), indicating that nonlinear effects in 
screening of external charges may be very strong. In particu- 
lar, this ratio reaches maximum values at certain distances r c 
that obviously depend on both the doping density |n| and the 
strength of external perturbation determined by Zq. [We note 
that the difference U- (r) — U+ (r) is always found to peak at 
the origin, r = 0.] 

The maxima in the ratios, seen in Fig. 3, may be explained 
by the fact that, for the hole doping (n < 0) of graphene in 
the presence of a positive external charge, there will be a local 
re-doping with electrons, or discharging of graphene, giving 
rise to a local shift of the tt electron band DOS, such that the 
condition U-(r c ) ~ HvFkF may be reached, indicating that 
the Fermi level is pushed back to cross the neutrality point at 
some distance r = r c . Since there are fewer states available 
in the DOS around the neutrality point, the screening ability 
of graphene is reduced around r = r c when n < 0, resulting 
in a higher value of the total potential than in the case of elec- 
tron doping (n > 0), so that one may expect that an inequality 
U-(r) > U+(r) > will hold for a range of distances r 
around r c . For example, in Fig. 3(a), the external charge is so 
close to graphene at zero temperature that it provides a strong 
enough perturbation, giving rise to the local discharging for 
all three doping densities, \n\ = 10 11 , 10 12 , and 10 13 cm~ 2 , 
so that three maxima in the ratio J7_ (r)/U + (r) occur around 
distances r c w 35.6, 12.7, and 4.8 A, respectively. The corre- 
sponding values of the potential f/_(r c ) at these distances are 
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FIG. 3: The ratio U-(r)/U+(r) of the nonlinear potential energies 
U- (r) and U+ (r) corresponding to, respectively, negative (hole dop- 
ing) and positive (electron doping) signs of the equilibrium charge 
carrier densities \n\ — 10 11 (solid [green] lines), 10 12 (dashed [red] 
lines), and 10 13 cm -2 (dash-dotted [blue] lines), is shown as a func- 
tion of the radial distance r (in A) in the plane of graphene for a 
proton at distances zo = 2 A with T — (panel a), Zo = 10 A with 
T = (panel b), and z = 10 A with T = 300 K (panel c), above 
graphene lying on an Si02 substrate with the gap h — 1 A. 



found to be 0.037, 0.137, and 0.495 eV, respectively, which 
scale reasonably close to the Fermi level shift at the three dop- 
ing densities, \e F \ = hv F k F w 0.037, 0.1 17, and 0.368 eV. 

On the other hand, when the charge is removed to distance 
zo = 10 A at zero temperature in Fig. 3(b), the perturba- 
tion is still strong enough to discharge graphene for the two 
lower doping densities [with the peaks occurring at similar 
distances, r c sa 47.8 and 15.0 A, and with similar potential 
values, U-{r c ) sa 0.041 and 0.153 eV, as in Fig. 3(a)], but is 
not sufficient do force the Fermi level to cross the neutrality 
point for the highest density of \n\ = 10 13 cm~ 2 , for which a 
maximal local discharging of graphene occurs directly under- 
neath the external charge. Furthermore, when the temperature 
is raised to T = 300 K for z = 10 A, the ratio f/_ (r) /U+ (r) 
for the two higher doping densities is barely affected, but the 
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FIG. 4: The relative error in the potential energy, U (r) = e<j>(r) (in 
%), from the nonlinear TF model for a proton at distance zo — 2 A 
from free, intrinsic (n — 0) graphene at zero temperature, due to the 
inclusion of the exchange and correlation energy 20 with values of the 
cutoff parameter being r\ = 1 (solid [red] line), 0.75 (dashed [green] 
line), and 0.5 (dotted [blue] line), as well as due to the nonlinear 
correction to graphene's n electron band density of state (dash-doted 
[black] line). 



ratio for the lowest density of |n| = 10 11 cm -2 appears to be 
largely suppressed in Fig. 3(c) as compared to Fig. 3(b). One 
can still see a maximum in this ratio around a distance similar 
to that in Fig. 3(b), i.e., r c w 37.2 A with U- » 0.045 eV, 
but the peak value of the ratio U-(r)/U+(r) for \n\ = 10 11 
cm" 2 has dropped from about 1.8 for T = K to about 1.2 
for T = 300 K. While the results in Fig. 3(c) confirm the 
conclusion drawn from Fig. 2 that, at room temperature, the 
screening ability of graphene is affected for sufficiently low 
doping densities, such that |n| < 10 11 cm -2 , it is now clear 
that the role of elevated temperature, when it prevails the ef- 
fects of doping density, is to the suppress the nonlinear effects. 

All results shown in Figs. 1-3 were obtained by taking into 
account in Eq. (Q]) the effects of nonlinearity in the band DOS 
of graphene, p(e), because we suspected that the value of 
the potential U(r) may exceed locally (that is, directly un- 
derneath the external charge) the cutoff value of about 1 eV 
that validates the linear approximation for p(e). Our calcu- 
lations show that the effect of this nonlinearity is relatively 
weak, giving corrections up to several percent for distances 
\zq\ > 1.5 A. This is illustrated in Fig. 4 for free, intrinsic 
(fi = 0) graphene at zero temperature with a charge Z = 1 
placed at zq = 2 A, where we show by the dash-dotted line the 
relative error in the total potential when Eq. (O is solved with 
density n from Eq. (O and from Eq. (Q]i with a nonlinear DOS 
p(s)£ One can see that the peak error of about 2 % occurs at 
the origin and diminishes at distances greater than a few A. 

We further estimate the effects of the exchange and correla- 
tion interactions, which have been neglected so far in solving 
the nonlinear TF equation (|9). We use the expression V xc (n) 
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FIG. 5: The effective dielectric constant e» in the image force, writ- 
ten as Fi m = (Ze/2zo) 2 (l/e* — 1), as a function of distance zo 
(in A) for a proton above free graphene at zero temperature. Re- 
sults from the nonlinear TF model are shown for equilibrium densi- 
ties n = (lower thick [black] solid line), ±10 12 (thin [red] solid 
and dashed lines, respectively), and ±10 13 cm" 2 (upper thick [blue] 
solid and dashed lines, respectively). Results from the linearized TF 
model and the RPA model are shown, respectively, by dash-dotted 
and dotted lines for densities \n\ — (lower thick [black] lines), 
10 12 (thin [red] lines), and 10 13 cm" 2 (upper thick [blue] lines). 



for the XC potential energy given by Polini at al. in the LDA 
and, since the formalism providing V xc (n) is restricted to in- 
trinsic graphene at zero temperature within the linear approx- 
imation for p(e)r& we solve the nonlinear equation, Eq. ( TTTb 
with Eq. ( fT2l . for free graphene (e£ g = 1) with the charge 

Z = 1 a distance zq = 2 A away. The result is compared 
to the solution when V xc is set to zero by showing in Fig. 4 
the relative error of such a comparison for several values of 
the cutoff parameter rj. 20 One can see in Fig. 4 that the rela- 
tive error due to the XC interactions is relatively small at short 
distances r, and is comparable to the error due to the nonlin- 
ear band DOS. However, the error due to the XC interactions 
increases and reaches a maximum of about 5% at distances on 
the order of r = 10 A or more, reverses its sign at still greater 
distances of about r = 100 A or more, and presumably con- 
tinues growing further in magnitude. While this is a relatively 
small error at radial distances where the total potential has a 
significant value, we note that the error due to the XC inter- 
action may be larger when external charge is placed further 
away from graphene, as noted earlier. 20 However, because of 
the limitation of the theory for XC interactions to local pertur- 
bations of charge carrier density relative to intrinsic graphene 
at zero temperature r^2°- we no longer pursue the analysis of 
the XC effects in our nonlinear TF approach. 

While the results in Figs. 1-4 elucidate local properties of 
the solution of the nonlinear TF equation, Eq. (0, we now 
turn to analyzing the image force F; m on a point charge as 
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a quantity that provides an integrated information on the ef- 
fects of doping and nonlinear screening in graphene. We first 
consider free graphene at zero temperature, and represent the 
nonlinear image force in the form reminiscent of the classi- 
cal image force of a point charge Ze in vacuum, a distance 
zq away from a layer of dielectric material with an effective 
dielectric constant e*, given by 



(Ze) 2 
4z 2 



1 



1 



(24) 



In this way, the z dependent parameter provides a mea- 
sure of the polarizability of free graphene. We use the same 
line styles and thicknesses as in Fig. 1 to show in Fig. 5 the 
results of the nonlinear TF calculations of as a function of 
zq for \n\ = 0, 10 12 , and 10 13 cm" 2 , along with the corre- 
sponding LTF and RPA results obtained from Eq. d22b with 
an appropriate polarization function by taking the derivative, 
Fi m = —dVim/dzQ. One can see in Fig. 5 a strong depen- 
dence of the nonlinear TF image force on both the magnitude 
and the sign of charge carrier density n, whereas the linear 
results seem to work only at large enough distances zq, with 
the RPA model showing a better agreement with the nonlin- 
ear TF results than the TF model. Notice that the slopes of 
the LTF lines follow from taking the derivative of the asymp- 
totic limit of the image potential in Eq. d23l l, and are given 
for n by the zero temperature limit of the inverse screen- 
ing length in Eq. (TE[ , q s — 4r s kp. On the other hand, the 
nearly horizontal lines for the nonlinear TF and the RPA mod- 
els with n — show that intrinsic graphene behaves as a layer 
of material with effective dielectric constants of rj 3.57 and 
« 1 + Trr s /2 ps 4.44, respectively. 

Finally, we analyze in Fig. 6 the image potential on a point 
charge Z = 1 above free graphene (panel a) and in the pres- 
ence of a Si02 substrate with zero gap (panel b), at zero tem- 
perature. We show the results due to the nonlinear TF and the 
RPA models for n — (thick solid and dotted lines, respec- 
tively) and ±10 13 cm" 2 (thin solid and dashed lines for the 
nonlinear TF, and thin dotted line for the RPA model), as well 
as the results due to the LTF model for \n\ = 10 13 cm -2 (thin 
dash-dotted line). We note that the nonlinear results were ob- 
tained by integrating the corresponding image force from zq 
up to typically 400 A. One notices a relatively close grouping 
of all results, indicating that the linear models provide good 
approximations, especially at high density and large distances 

Zq. 

However, the effects of doping of graphene are seen to 
be still quite strong giving, e.g., in the nonlinear TF model 
for free graphene the image potential of Vj m w — 0.32 eV 



10 A when 



10 cm , as opposed to 



at Zq 

Vi m ~ —0.26 eV found at the same distance above intrin- 
sic graphene. This points to possibly strong effects of doping 
in the asymptotic region of distances of relevance to the im- 
age potential states. 28 While the discrepancy between the RPA 
and the nonlinear TF results, seen in Fig. 6 for free graphene 
at zero doping, stems from the difference seen in Fig. 5 be- 
tween the effective dielectric constants of intrinsic graphene 
in those two models, one notices a near-perfect agreement 
of the RPA model with the nonlinear TF model in graphene 
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FIG. 6: Image potential Vim (in eV) of a proton as a function of 
its distance zo (in A) above graphene at zero temperature, for free 
graphene (upper panel) and for graphene lying on an Si02 substrate 
with zero gap (lower panel). Results from the nonlinear TF model 
are shown for equilibrium densities n — (upper thick [black] solid 
line) and ±10 13 cm" 2 (lower thin [blue] solid and dashed lines, re- 
spectively). Results from the RPA model are shown by dotted lines 
for densities |n| = (upper thick [black] line) and 10 13 cm -2 (lower 
thin [blue] line), as well as from the linearized TF model for density 
Inl = 10 13 cm -2 (thin [blue] dash-dotted line). 



doped by electrons to n = 10 13 cm" 2 . However, nonlin- 
ear effects are still quite strong, especially at short distances, 
as illustrated by the observed asymmetry in the nonlinear TF 
model with respect to the sign of n ^ 0. For example, one 
finds in Fig. 6(a) that the image potential takes the value of 
Vim ~ — 1.93 eV at zo ~ 1.5 A above free graphene with 
n = 10 13 cm" 2 , as opposed to V5 m « —1.66 eV at the same 
distance with n — — 10 13 cm -2 . This asymmetry due to dop- 
ing of graphene by electrons or holes may have interesting 
and important consequences for, e.g., chemisorption of a Li 
atom, where the image potential shift of its 2s orbital level 
may be controlled by the applied gate potential and used to 
move around the resonance in the local DOS, and even possi- 
bly break the ionic bond between the Li atom and graphene. 
Finally, we note that we have estimated numerically the ef- 
fects of non-zero temperature and the XC interactions in the 
nonlinear image potential for intrinsic graphene, and found 
that both these effects are negligible compared to the above 
effects of the doping density and nonlinear screening. 
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IV. CONCLUDING REMARKS 

We have solved a nonlinear TF equation for the radial de- 
pendence of electric potential in the plane of single-layer 
graphene due to an external point charge in the presence of 
a dielectric substrate with a finite graphene-substrate gap, h, 
paying special attention to the effects of equilibrium charge 
carrier density n, temperature T, and separation between the 
charge and graphene \zq\. Large effects were found due to 
variations in both the magnitude and the sign of n, illustrating 
the importance of both doping of graphene and the nonlin- 
ear screening, respectively. Temperature was found to mostly 
affect screening at low doping densities, satisfying the in- 
equality kp = y/Tr\n\ < ksT/ (Hvf), in such a way as to 
suppress the nonlinear effects. In addition, the existence of 
a non-zero gap, h, between the substrate and graphene was 
found to exert non-negligible effects on the potential, mostly 
at short radial distances. We have moreover analyzed the ef- 
fects in the potential due to nonlinear corrections in the den- 
sity of states of graphene's ir electron bands, as well as due 
to the exchange and correlation interactions for the case of 
free, intrinsic graphene at zero temperature. While the for- 
mer effect gives corrections of up to a few percent at positions 
directly underneath the external charge and diminishes at dis- 
tances further out, the latter effect gives rise to the corrections 
of up to 5 % at intermediate and large radial distances. 

Comparisons were made with the results from a linearized 
TF (LTF) equation and from the RPA model of dielectric 
screening in graphene. While the LTF results are generally 
close to the nonlinear TF results at large radial distances and 
high densities \n\ only, the RPA model also exhibits an im- 
proved agreement with the nonlinear TF model at short radial 
distances, owing to the short wavelength dielectric constant 
of graphene, which results from the inter-band electron transi- 
tions captured by the RPA modeli22i2& Unlike the TF models, 
the RPA results exhibit Friedel oscillations around the poten- 
tial from the linearized TF model at large radial distances in 
doped graphene, with amplitudes that increase with increas- 
ing gap h, but are dampened by increasing separation \ zq \ and 
increasing temperature. 

Our most important conclusion is that nonlinear effects are 
strong over a broad range of radial distances, even at high 
doping densities |n| and large separations \zq\, as illustrated 
by the large ratios of the potential evaluated from the non- 
linear TF model with the same amounts of doping by holes 
(n < 0) and by electrons (n > 0). This may be explained by 
a local shift of graphene's density of states so that the Fermi 
level is forced to cross the neutrality point in that density at 
a certain radial distance, thereby reducing graphene's polariz- 



ability when doping occurs with carriers of the same charge 
sign as the external particle. This asymmetry in the scatter- 
ing potential for charge carriers in graphene with respect to 
the sign of n may be responsible for the observed asymme- 
try in graphene's conductivity as the sign of the gate poten- 
tial changes. 46 However, such an effect of nonlinear screening 
of external charges will be suppressed at low doping densi- 
ties when the temperature is sufficiently elevated, as described 
above. 

Finally, we have analyzed the image interaction of an ex- 
ternal charge due to polarization of graphene, where we com- 
pared the results evaluated from the solution of the nonlinear 
TF equation with those from the LTF and the RPA models. 
After elucidating the strong doping and nonlinear effects in 
the image force above free graphene at zero temperature, we 
have presented results for an image potential obtained by nu- 
merical integration of the nonlinear image force up to large 
distances from graphene, and compared them with the results 
of the linear models. The nonlinear image potential was found 
to exhibit relative variations due to doping of graphene up to 
\n\ = 10 13 cm~ 2 , which may reach over 20 % at distances 
\zq\ ~ 10 A, as well as due to the nonlinear screening, where 
relative variation with the sign of n may reach about 20 % at 
short distances, on the order of \zq\ ~ 1 A. These variations 
in the image potential were found to be somewhat reduced in 
the presence of an Si02 substrate. 

Our results for the electric potential in the plane of graphene 
due to an external charge may be relevant for calculations of 
its conductivity based on the Boltzmann transport modelriS 
where this potential may be used directly in an expression 
for the transport relaxation time in the Born approximation, 
to reveal the effects of doping, nonlinear screening and tem- 
perature on conductivity. While this task is left for a future 
contribution, we comment here that our nonlinear TF results 
are likely to yield calculable effects due to the asymmetry in 
charge of the external particles, 46 based on the presently ob- 
served asymmetry with respect to the sign of n for a posi- 
tive external charge. Moreover, our results for the nonlinear 
image potential may be found helpful in studying chemical 
processes near graphene, e.g., alkali atom chemisorption and 
intercalation, 24 as well as in the recent work on the electron 
image-potential states near graphene^ 
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